/*
  Copyright (C) 2011, Kyungjoo Kim
  
  This file is part of LINAL (LINear ALgebra)
  
  All rights reserved.

  This library is free software; you can redistribute it and/or
  modify it under the terms of the GNU Lesser General Public
  License as published by the Free Software Foundation; either
  version 2.1 of the License, or (at your option) any later version.
  
  This library is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  Lesser General Public License for more details.
  
  You should have received a copy of the GNU Lesser General Public
  License along with this library; if not, write to the Free Software
  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
  Also add information on how to contact you by electronic and paper mail.
  
  
  Kyungjoo Kim
  iamkyungjoo@gmail.com
*/
#include "gpu_perform.hxx"

int main(int argc, char **argv) {

  if (argc != 7) {
    printf("Try :: chol [n_thread] [n_gpu] [blocksize] [cpu_blocksize] [n_dof] [n_itr]\n");
    return 0;
  }

  // ---------------------------------------
  // ** Initialization
  int n_dof, blocksize, n_thread, n_itr, n_gpu, cpu_mn, n_block;
  n_thread  = atoi( (argv[1]) );
  n_gpu     = atoi( (argv[2]) );
  blocksize = atoi( (argv[3]) );
  cpu_mn    = atoi( (argv[4]) );
  n_dof     = atoi( (argv[5]) );
  n_itr     = atoi( (argv[6]) );
  n_block   = 64;

  printf("** TEST ENVIRONMENT **\n");
  printf("NDOF      = %d\n", n_dof);
  printf("Blocksize = %d\n", blocksize);
  printf("N thread  = %d\n", n_thread);
  printf("N gpu     = %d\n", n_gpu);
  printf("CPU block = %d\n", cpu_mn);
  printf("Iteration = %d\n", n_itr);

  int b_mn[2];
  b_mn[0] = b_mn[1] = blocksize;

  FLA_Init();

  double flop = linal::get_flop_chol( 0, n_dof );

  // ---------------------------------------
  // ** Matrices
  double 
    t_base, t_flash_repack, t_linal_repack, t_temp, 
    t_flash_decompose, t_linal_decompose;

  linal::Flat_ A;
  A.create(TEST_DATATYPE, n_dof, n_dof);
  FLA_Random_matrix( ~A );

  FLA_Obj hA_fla;
  t_base = FLA_Clock();
  FLASH_Obj_create_hier_copy_of_flat(~A, 1, (dim_t*)b_mn, &hA_fla);
  t_flash_repack = FLA_Clock()-t_base;

  linal::Hier_ hA_linal;
  t_base = FLA_Clock();
  hA_linal.create(A, blocksize, blocksize);
  t_linal_repack = FLA_Clock()-t_base;
  
  // ---------------------------------------
  // ** FLASH
  FLASH_Queue_set_num_threads(n_thread);
  FLASH_Queue_set_sorting(TRUE);
  FLASH_Queue_set_caching(TRUE);
  FLASH_Queue_enable_gpu();
  FLASH_Queue_set_gpu_num_blocks(n_block);

  t_flash_decompose = MAX_TIME;
  for (int q=0;q<n_itr;++q) {
    printf("*** FLASH::CHOL BEGIN ***\n");
    t_base = FLA_Clock();

    FLASH_Queue_begin();
    FLASH_Chol(FLA_LOWER_TRIANGULAR, hA_fla);
    FLASH_Queue_end();

    t_temp = FLA_Clock()-t_base;
    printf("*** FLASH::CHOL END ***\n");

    t_flash_decompose = min(t_temp, t_flash_decompose);
  }

  FLASH_Queue_disable_gpu();

  // ---------------------------------------
  // ** LINAL
  switch (TEST_COMPUTING) {
  case LINAL_CPU_GPU:
    omp_set_num_threads(n_thread+n_gpu);
    break;
  default:
    omp_set_num_threads(n_thread);
    break;
  }
  int n_item_bin = n_block, n_item_work = 1;
  linal::init_gpu(n_thread, n_gpu, n_item_bin, n_item_work);
  linal::set_computing_model( TEST_COMPUTING );
  linal::set_cpu_blocksize( cpu_mn );
  linal::create_flat_gpu( TEST_DATATYPE, b_mn[0], b_mn[1]);

  t_linal_decompose = MAX_TIME;
  for (int q=0;q<n_itr;++q) {
    printf("*** LINAL::CHOL BEGIN ***\n");
    t_base = FLA_Clock();
#pragma omp parallel 
    {
#pragma omp single nowait
      linal::dense::chol(FLA_LOWER_TRIANGULAR, hA_linal);
    }
    t_temp = FLA_Clock()-t_base;
    printf("*** LINAL::CHOL END ***\n");

    t_linal_decompose = min(t_temp, t_linal_decompose);
  }

  int cache_hit, gpu_hit;
  linal::get_gpu_data_access_counter(cache_hit, gpu_hit);
  printf("cache hit %d, total gpu hit %d, ratio %lf\n",
         cache_hit, gpu_hit, ((double)cache_hit)/((double)gpu_hit));

  linal::set_computing_model( LINAL_CPU );
  linal::free_flat_gpu();
  linal::finalize_gpu();

  switch (TEST_COMPUTING) {
  case LINAL_CPU_GPU:  printf("CPU + GPU BLAS\n");    break;
  case LINAL_GPU:      printf("GPU BLAS\n");          break;
  default:	       printf("CPU ONLY\n");          break;
  }
	
  printf("----------------------------------------------\n");
  printf("*** Report Chol ***\n");
  printf("N dof       = %d\n", n_dof);
  printf("Blocksize   = %d\n", blocksize);
  printf("N thread    = %d\n", n_thread);
  printf("N GPU       = %d\n", n_gpu);
  printf("N iteration = %d\n", n_itr);
  printf("----------------------------------------------\n");
  printf("Time Supermatrix   = %6.3lf [sec]\n", t_flash_decompose);
  printf("Time LINAL         = %6.3lf [sec]\n", t_linal_decompose);
  printf("----------------------------------------------\n");
  printf("FLOPS Supermatrix  = %6.3lf [Gflops]\n", 
	 flop/t_flash_decompose/1.0e9);
  printf("FLOPS LINAL        = %6.3lf [Gflops]\n", 
	 flop/t_linal_decompose/1.0e9);
  printf("----------------------------------------------\n");
  printf("With repacking cost \n");
  printf("FLOPS Supermatrix  = %6.3lf [Gflops]\n", 
	 flop/(t_flash_decompose+t_flash_repack)/1.0e9);
  printf("FLOPS LINAL        = %6.3lf [Gflops]\n", 
	 flop/(t_linal_decompose+t_linal_repack)/1.0e9);
  printf("----------------------------------------------\n");
  
  // ---------------------------------------
  // ** Matrix
  FLASH_Obj_free(&hA_fla);
  hA_linal.free();
  A.free();

  printf("*** TEST FINISHED ***\n");

  // ---------------------------------------
  // ** Finalization
  FLA_Finalize();
  return 0;
}
